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This  section  contains  a  factual  listing  of  publications,  research  lectures,  and  personnel 
supported  in  whole  or  in  part  under  ONR  Grant  No.  N00014-90-J-4071.  The  narrative  to 
follow  in  Section  II  will  be  keyed  to  the  item  numbers  in  these  lists.  Copies  of  the  abstracts  of 
published  and  submitted  papers  will  follow  in  Section  III. 


A.  Papers  Published  in  Refereed  Journals 

1.  Overfell,  P.L.  and  C.S.  Kenney,  "Comparison  of  the  Propagation 
Characteristics  of  Bessel,  Bessel-Gauss,  and  Gaussian  Beams  Diffracted  by  a 
Circular  Aperture,"  J.  Optical  Soc.  of  America  A,  Vol.  8,  No.  5,  May  1991. 

2.  Kenney,  C.S.  and  P.L.  Overfelt,  "A  Simple  Approach  to  Mode  Analysis  for 
Parabolic  Waveguides,"  IEEE  Trans.  Microwave  Th.  and  Tech.,  Vol.  39,  No. 
3,  March  1991. 

3.  Hench,  J.J.,  C.S.  Keruiey,  and  A.J.  Laub,  "Methods  for  the  Numerical 
Integration  of  Hamiltonian  Systems,"  to  appear  in  Circuits,  Systems,  and 
Signal  Processing,  1994. 
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F.  Contributed  Presentations  at  Scientific/Technical  Society  Conferences 

1.  Helical  Solutions  to  Fiber-Payout  Problems:  for  the  Fifth  Tri-Service 
Conference  on  Fiber  Optic  Guided  Missiles,  Naval  Ocean  Systems  Center,  San 
Diego,  CA,  November  5-7,  1990. 

2.  Inner  and  Outer  Modes  for  Fiber  Dynamics:  for  the  Sixth  Tri-Service 
Conference  on  Fiber  Optic  Guided  Missiles,  Naval  Postgraduate  School, 
Monterey,  CA,  March  18-21,  1991. 


G.  Technical  Reports  Published  or  Non-refereed  Journals 

1.  Hench,  J.J.,  C.S.  Kenney,  and  A.J.  Laub,  "A  Pad6-Magnus  Method  for  the 
State  Transition  Calculation  of  a  Linear  Time-Varying  System,"  Technical 
Report  SCL  91-03,  University  of  California,  Santa  Barbara,  Electrical  and 
Computer  Engineering  Dept.,  ^ptember,  1991. 

2.  Kenney,  C.S.,  "Inner  and  Outer  Wave  Modes  For  Fiber  Dynamics,"  Naval 
Weapons  Center  Technical  Paper  No.  NWC  TP  7149,  June  1991. 


H.  Patents 

1.  Hoban,  F.,  G.  Hewer,  D.  Harms,  and  C.  Kenney,  patent  applied  for  "Internal 
Fiber  Optics  Payout  System  with  Rotary  Drive  and  Control,"  Navy  Case  No. 
73679,  June  5,  1991. 
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n.  Narrative 


In  this  section,  we  provide  a  summary  of  the  research  conducted  during  the  period 
September  1990  to  September  1993.  Because  of  an  administrative  oversight,  ftmding  for  this 
grant  covered  only  the  period  September  1990  through  September  1991.  All  research 
conducted  subsequent  to  September  1991  was  unfunded.  A  more  detailed  description  of  some 
of  our  research  topics  is  contained  in  the  original  proposal  for  the  grant. 
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The  main  goal  of  our  work  was  to  develop  sound  numerical  algorithms  for  solving  the 
partial  differential  equations  that  arise  in  string  dynamics  problems,  especially  of  the  type 
encountered  by  fiber  optic  guided  missiles.  In  pursuit  of  this  goal,  we  worked  closely  with 
engineers  (Fay  Hoban  and  others)  and  mathematicians  (Gary  Hewer)  and  physicists  (Pam 
Overfelt)  involved  in  the  SKYRAY  project  at  the  Naval  Weapons  Center  in  China  Lake, 
California.  SKYRAY  was  a  demonstration  project  for  fiber  payout  at  high  velocities  (8(X) 
ft/sec  and  greater)  and  this  cooperative  effort  resulted  in  two  patent  applications;  a  novel  inside 
payout  design  and  a  rotating  outer  payout  design  (see  H.  1  and  H.2'  respectively  above  and  the 
discussion  below). 

The  partial  differential  equations  that  describe  fiber  payout  (see  the  grant  proposal)  reduce 
to  ordinary  differential  equations  subject  to  boundary  conditions  in  the  case  of  uniform  rotary 
motion.  This  reduction  served  as  the  point  of  departure  for  our  study  for  several  reasons. 
First,  the  assumption  of  uniform  rotary  motion  is  very  close  to  that  measured  experimentally, 
and  at  the  same  time  reduces  the  number  of  system  parameters  to  a  small  group  that  can  be 
easily  varied  for  design  studies.  Secondly,  the  resulting  ODE  BVP  can  be  solved  exactly  in 
terms  of  elliptic  sine  ftmctions  for  the  special  case  of  inner  payout  with  zero  air  resistance  (see 
G.2  above).  This  is  significant  because  the  inner  payout  modes  are  very  insensitive  to  changes 
in  air  resistance,  and  consequently  can  be  well  modeled  by  the  zero  resistance  case.  Moreover, 
the  inner  payout  boundary  conditions  can  be  used  (in  a  half  wave  number  mode)  to  generate 
approximate  boundary  conditions  for  the  outer  payout  modes.  The  vector  form  of  Newton’s 
method  then  provides  a  way  of  finding  the  true  boundary  conditions  for  the  outer  modes.  The 
details  describing  this  procedure  were  developed  in  a  series  of  eight  unpublished  working  notes 
for  the  design  group  at  NWC  and  summarized  in  G.2.  Subsequently  we  implemented  these 
procedures  as  MATTLAB  macros  in  a  software  package  called  FIBERPACK.  During  the 
testing  of  this  package  we  discovered  that  the  inner  payout  modes  produced  far  less  tension  in 
the  fiber  compared  to  the  outer  modes  especially  during  high-altitude  simulations.  We  also 
found  that  rotating  the  bobbin  to  assist  payout  reduces  the  tension  even  further.  Reductions  of 
this  kitxl  are  crucial  to  high-speed  payout  because  the  tension  increases  as  the  square  of  the 
payout  velocity.  These  observations  formed  the  basis  of  the  rotating  and  inner  payout  design 
patent  applications  mentioned  above. 

Problems  related  to  signal  attenuation  during  fiber  transmission  were  also  considered. 
Working  with  Pam  Overfelt  of  the  Physics  Division  at  NWC,  we  developed  a  variety  of 
methods  for  analyzing  and  optimizing  power-handling  characteristics  as  well  as  modes  of 
propagation  (see  A.l,  A.2,  and  Cl).  The  principal  benefit  of  these  studies  has  been  an 
increased  awareness  of  the  signal  attenuation  caused  by  sharp  kinks  or  bends  in  the  fiber. 
Such  kinks  appear  to  be  exaggerated  by  torsion  effects  related  to  payout  designs  which  use 
constraining  shrouds  as  described  below. 

In  future  research  it  would  be  desirable  to  extend  our  effort  to  include  problems  in  which 
the  fiber  contacts  a  constraining  shroud  and  may  suffer  time-varying  impluse  forces.  The 
shroud  problem  is  of  interest  because  of  the  need  to  force  the  fiber  in  certain  directions,  such 
as  away  from  the  rocket  plume.  The  shroud-fiber  contact  point  induces  a  torsion  effect  into  the 
governing  differential  equations  that  makes  numerical  computation  difficult  and  is  sure  to 
produce  some  interesting  results. 

Time-varying  impulse  forces  also  need  to  be  studied  because  the  payout  environment  is 
subject  to  extreme  aerodynamic  forces  and  because  the  design  engineers  at  NWC  want  to  use 
high  frequency  dither  to  reduce  frictional  effects  at  the  shroud  contact  point.  The  simplest 
models  of  these  forces  result  in  ordinary  differential  equations  with  time-varying  coefficients. 
Fortunately  we  were  able  to  complete  an  initial  study  of  this  problem  in  A.3  and  it  appears  that 
very  efficient  solution  procedures  are  possible. 


-  5  - 

We  are  excited  about  the  progress  that  was  made,  but  many  important  open  research 
problems  remain. 

m.  Abstracts  of  Published  and  Submitted  Papers 
See  attached. 
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W«  use  the  scalar  Kitchhoff-Huygens  diffraction  integral  to  obtain  analytic  expreuiont  for  both  axial  and  trans¬ 
verse  intensity  distributions,  assuming  normal  incidence  on  a  circular  aperture  for  four  types  of  incident  field;  1 1 ) 
plane  wave.  (2)  Bessel  beam.  (3)  Gaussian  beam,  and  14)  Bessel-Gauss  beam.  We  use  the  Fresnel  approximation  to 
obtain  the  axial  intensity  as  a  function  of  distance  from  the  aperture.  We  consider  both  Fresnel  and  Fraunhofer 
diffraction  for  the  case  of  the  transverse  intensity  distributions.  For  the  axial  case,  we  And  that  the  Bessel-Gauss 
beam  performs  worse  than  the  Bessel  beam,  in  terms  both  of  the  magnitude  of  intensity  and  of  its  ability  to  extend  a 
distance  from  the  aperture.  In  the  transverse  cose,  we  find  that  the  Bessel-Gauss  beam  performance  in  terms  of 
remaining  nearly  diffraction  free  over  a  given  distance  is  highly  dependent  on  the  relationship  among  the  aperture 
radius,  the  beam  waist  parameter,  and  the  transverse  wave  number. 


1.  INTRODUCTION 

Recently,  there  haa  been  much  interest  in  the  possible  prop¬ 
agation  of  diffraction-free  beams  at  optical  frequencies.'*' 
These  include  the  Bessel  beam'  *  and  the  Bessel-Gauss 
beam. '  In  previous  work  the  Bessel-Gauss  beam  has  been 
considered  only  for  free-space  propagation,  while  the  Bessel 
beam  was  analyzed  by  numerical  simulation  and  verified 
subsequently  through  experiment. 

In  what  follows  we  use  scalar  Kirchhoff  theory  to  obtain 
analytic  expressions  for  the  axial  and  transverse  intensity 
distributions,  assuming  normal  incidence  upon  a  circular 
aperture  for  four  types  of  Held:  (1)  plane  wave,  (2)  Bessel 
beam,  (3)  Gaussian  beam,  and  (4)  Bessel-Gauss  beam. 
Since  the  Bessel-Gausa  beam  is  the  most  general  of  the 
incident  fields  considered  above,  once  ita  intensity  distribu¬ 
tion  is  found  the  distributions  of  the  remaining  fields  can  be 
recovered  by  using  limiting  cases  of  the  Bessel-Gauss  inten¬ 
sity  function.  The  analysis  for  the  axial  intensity  is  shown 
in  Section  2,  and  that  for  the  transverse  intensity  is  shown  in 
Section  3.  In  Section  4  we  discuss  the  Fraunhofer  diffrac¬ 
tion  pattern  generated  by  an  incident  Bessel  beam  and  show 
that  the  amount  of  energy  conteined  within  any  specified 
radius  can  be  controlled  by  changing  the  value  of  the  trans¬ 
verse  wave  number,  which  governs  the  sharpness  of  the  cen¬ 
tral  maximum  of  the  Bessel  beam.  Section  d  contains  our 
conclusions. 

2.  AXIAL  INTENSITY  IN  THE  FRESNEL 
REGION 

Assuming  exp(iut)  time  dependence  throughout,  cylindrical 
coordinates  (p,  a.  r)  and  incident  fields  that  are  cylindrically 
symmetric  (independent  of  o).  the  scalar  Kirchhoff-Huy- 
gens  diffraction  integral  (for  normal  incidence)  is  given  by* 


iko  exp(-i*o«) 


exp(-ikiyR) 


(U 


where  A  is  the  aperture  in  a  plane  screen,  *  2jr/\,,  is  the 
free-space  wave  number,  ^inr(p'>  0)  is  the  incident  field  at  c 
■  0,  and 

R  ■  (p'^  +  p*  +  z*  -2pp'  cos(0'  -  <&)]'  -.  (2) 

For  the  case  of  axial  intensity  as  a  function  of  distance  from 
the  aperture  in  the  upper  half-plane,  r  >  0  (see  Fig.  1 ).  we  set 
p  w  0.  Thus  Eq.  (2)  becomes 

R  -  (p'*’ +  13) 


R  •  z  +  ^, 

2z 


14) 

15) 


by  using  the  Fresnel  approximation  for  R;  i.e..  R  «  z  in  the 
denominator  and  R  given  by  relation  (5)  in  the  exponential 
portion  of  Eq.  (1).  The  scalar  diffracted  field  given  by  Eq. 
(1)  becomes 

(k„  exp(~ik,yi} 

*(0.  z)  -  - - 

2ir  z 

xi  I  \t,„,(p'.  0)exp(-iA„p'72z)p'dp'd4)  i6) 

in  .n 

for  a  circular  aperture  of  radius  a  in  a  plane  screen.  .At  this 
point,  we  should  note  that  the  Raleigh-Sommerfeld  integral 
(in  the  paraxial  limit)  would  have  led  to  Eq.  16)  also  las  was 
pointed  out  by  an  unknown  reviewer).  We  want  to  evaluate 
Eq.  (6)  for  four  incident  fields:  plane  wave.  Bessel  beam. 
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A  Simple  Approach  to  Mode  Analysis  for 
Parabolic  Waveguides 

Charles  S.  Kenney,  Associate  Member,  IEEE,  and  P.  L,  Overfelt,  Member.  IEEE 


Abstract  — Difllcuity  in  obtaining  accurate  values  for  parabolic  crlin- 
der  hinctions  has  been  an  impediment  to  mode  analysis  for  parabolic 
waveguides.  A  simple  method,  based  on  one^imensional  analytic  con¬ 
tinuation.  is  presented  which  gives  essentiaily  exact  values  for  these 
functions:  i.e..  the  relative  error  in  the  computed  result  is  on  the  order  of 
the  machine  round>oir.  When  supplemented  with  a  Newton-Poisson 
shooting  method  and  simple  homotopy  techniques,  this  continuation 
method  can  he  used  to  find  the  TE  and  TM  mode  eigenvalues,  and 
associated  separation  constants,  for  arbitrary  parabolic  domains.  These 
methods  arc  then  used  to  compute  a  power  handling  efficiency  factor  for 
a  range  of  parabolic  regions. 

I.  Introduction 

WE  consider  parabolic  cylinders  of  uniform  cross  section 
in  the  confocal  parabolic  coordinates  (f.n.;).  which 
are  related  to  rectangular  coordinates  (X.Y.Z)  via 

.»Z  (I) 

(see  Fig.  1).  The  cross  sections  of  interest  consist  of  the 
interior  regions  ft  *  ftffi.rjn)  bounded  on  the  right  by  the 
curve  Tj  -  Tj„  and  on  the  left  by  the  curve  f  *  f„. 

Assuming  a  uniform,  perfectly  conducting  waveguide  of 
parabolic  cross  section  with  c dependence,  we  use 


where  ilt  satisfies 

\l/x.t 'Ityy  “  -  k-\b  in  ft  (2) 


subject  to  the  boundary  conditions. 

£, -0 

on  lift  (TM  modes) 

(3) 

or 

1 

o 

on  sn  (TE  modes). 

(■») 

Here  (^ft  denotes  the  boundary  of  ft.  ri/fln  is  the  outward 
normal  derivative,  and  variable  subscripts  indicate  different!' 
ation.  In  <2).  k-  *  kr,  -  13'  and  k,]  -  ai-e„fs„. 
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Expressing  (2)  in  parabolic  coordinates  gives 

<lt((  +  -  k‘{^^ V')'!'-  (5) 

Using  separation  of  variables  with  \li  «  U’lf  )1^(  rj)  we  tind 

Uf(^(k'f'-^a)U~0  (6) 

|/,,+(*-r7--Qi)l/-0  (') 

where  a  is  the  separation  constant.  Equations  (6)  and  (') 
must  be  supplemented  by  the  boundary  conditions 

£/(fii)-0  l^(T7n)-0  (TM  modes)  (X) 


i>{(fii)“0-  K,(n())*0  (TE  modes)  (9) 

Although  (6)  and  (7)  can  be  solved  for  any  a  and  k.  the 
boundary  conditions  (8)  and  (9)  are  satisfied  only  at  discrete 
pairs  ia.k),  when  („  and  are  fixed. 

In  the  following  section,  we  show  that  solutions  to  (6)  and 
(7)  can  be  computed  via  one-dimensional  analytic  continua¬ 
tion.  Section  III  discusses  a  Newton-Poisson  shooting 
method  for  finding  the  separation  constants,  a.  and  eigenval¬ 
ues.  k.  for  fixed  f„  and  tj„.  This  method  is  easy  to  use  and 
reveals  some  inaccuracies  in  previously  published  work.  In 
particular.  Tables  I  and  II  give  pairs  of  values  (a.t))  «ia  /  2k. 
v/2r )  to  seven  significant  digits  for  the  case  f ,  -  n.i.  and 
show  that  some  entries  in  similar  tables  from  [IJ  have  onK 
one  digit  of  accuracy.  The  values  in  Tables  I  and  II  were 
used  to  generate  (via  a  simple  homotopy  method)  the  eigen- 
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Methods  for  the  Numerical  Integration 
of  Hamiltonian  Systems* 

J.J.  Hench,  C.S.  Kenney,  and  A.J.  Laub 
Department  of  Electrical  and  Computer  Engineering 
University  of  California 
Santa  Barbara,  CA  93106-9560 


Abstract 

In  order  to  compute  an  infinite  horizon  optimal  controller  for  a  linear  periodic  system  via  an 
invariant  subspace  method,  the  computation  of  the  period  map  associated  with  the  Hamilton- 
Jacobi- Bellman  equations  is  required.  In  this  paper  we  discuss  methods  for  the  numerical 
integration  of  such  Hamiltonian  systems.  Two  numerical  integration  techniques  are  introduced. 
A  method  is  developed  whereby  symplectic  invariants  associated  with  the  Hamilton-Jacobi- 
Bellman  equations  are  preserved.  Also,  a  shifting  scheme  is  introduced  which  in  effect  swaps 
the  roles  of  the  stable  and  unstable  invariant  subspaces  by  using  the  semigroup  property  of 
state  transition  matrices.  A  shift  is  introduced  into  the  resultant  initial  value  problem  assuring 
that  the  eigenvalues  of  a  differential  equation  reside  in  the  region  of  absolute  stability  for  an 
appropriate  numerical  integration  routine.  These  techniques  are  then  compared  to  standard 
numerical  integration  routines  to  ascertain  their  efficiency  and  accuracy. 

Keywords:  Fade  approximation,  linear  systems,  periodic  systems,  periodic  Riccati  equation, 
state  transition  matrix,  invariant  subspaces. 

AMS(MOS)  subject  classification:  65L05 

Abbreviated  title:  Numerical  Integration  of  Hamiltonian  Systems. 


1  Introduction 


In  this  paper,  we  examine  various  methods  of  computing  the  optimal  feedback  controller  gains  for 
linear,  periodic  time- varying  systems 

i(0  =  A{t)x{t)  +  B{tMt)  (1) 

with  quadratic  cost  functional: 

mm  ^  J  -I-  u^/Z(t)u|  dt  (2) 

‘This  research  was  supported  by  the  Air  Force  Office  of  Scientific  Research  under  Grant  No.  AFOSR-91-0240, 
the  National  Science  Foundation  under  Grant  No.  ECS-9120643,  and  the  Office  of  Naval  Research  under  Grant  No. 
N00014-90-J-4071  and  No.  N00014-92-J-1706. 
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Optimizing  waveguide  ctoss-  sections 
with  respect  to  power  handling  capability 


C.S.  Kennej^.  P.L.  Overfelt* 

1.  Department  of  Electrical  and  Computer  Engineering.  University  of  California. 

Santa  Barbara.  CA-93106  2.  Pities  Division.  Research  Department 
Naval  Weapons  Centre.  China  Lake.  CA*9355S*6001.USA. 


ABSTRACT 

Ao  »n«lysia  ta  pvaa  of  the  power-headlios  ceps- 
bilitiee  of  four  familiei  of  waveguide  eroaa'teetiooal 
shapes,  including  recUngular,  elliptical,  parabolic 
and  rhombic.  Using  an  efficiency  factor  based  on 
the  two  lowest  TE  modes,  it  is  seen  that  ellipses 
are  slightly  better  than  the  standard  2:1  rectangu¬ 
lar  waveguide  and  that  for  the  waveguides  consid¬ 
ered  here,  power-  handling  is  maximised  when  the 
second  and  third  TE  mode  eigenvalues  are  eqnsL 
Since  this  efficiency  factor  is  relatively  insensitive  to 
slight  variations  in  cross-sectional  shape,  good  esU- 
mates  of  power-handling  capability  can  be  made  by 
comparison  with  rectangular  shapes  or  other  simple 
geometries  of  the  same  aspect  ratio. 

INTRODUCTION 

This  paper  addresses  the  problem  of  finding  tbs 
waveguide  crosa-seetional  shape,  (1,  which  optimises 
power  handUng  capability  as  measured  by  the  effi¬ 
ciency  factor,  7,  defined  by  Baum  (1988).  This  gen¬ 
eral  problem  is  difficult  because  7  depends  on  the 
first  two  TE  modes  of  (1,  i.e.  the  first  tvro  eigenval¬ 
ues  of  the  Laplacian  over  (1  with  Neumann  bound¬ 
ary  conditions: 

Atlr «  in  0.  (1) 


>0  on  Ml,  (2) 

where  Ml  denotes  the  boundary  of  (1,  d/Mi  is  the 
outward  normal  derivative  and  0  <  hf  <  h)  <  ■  ■ 
More  specifleally,  the  efficiency  fhetor  7  »  7<n)  is 
deffned  as 


Since  ansJytieal  espressioas  for  hi,  hr  and 
are  not  known  for  gsasral  domains,  we  restrict  our 
attention  to  four  families  of  shapes  (rectangular, 
elliptical,  parabolic,  and  rhombic)  which  can  be 
parametria^  by  a  single  variable,  r.  For  each  of 
these  shapes,  r  is  the  ratio  of  the  height  to  the 
width,  which  we  refer  to  as  the  aspect  ratio.  For 
recungular,  elliptical  and  parabolic  regions,  exact 
values  of  7  w  7(r)  were  found  by  using  respec¬ 
tively,  trigonometric  functions,  Mathien  fimetions 
and  parabolic  cylinder  functions.  For  rhombic  re¬ 
gions,  approximate  efficiency  factors  were  deter¬ 
mined  by  using  highly  accurate  collocation  methods 
u  described  by  OveAlt  and  Kenney  (1990). 

There  are  several  notable  features  which  can  be 
drawn  from  this  study.  First,  for  the  four  types 
of  regions  considered  here,  the  efficiency  factor  is 
15  maximixed  when  the  second  and  third  eigenvalues 


